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We present a study of the isotropic-nematic interface in a system of hard spherocylinders. First 
we compare results from Monte Carlo simulations and Onsager density functional theory for the 
interfacial profiles of the orientational order parameter and the density. Those interfacial properties 
that are not affected by capillary waves are in good agreement, despite the fact that Onsager theory 
overestimates the coexistence densities. Then we show results of a Monte Carlo study of the capillary 
waves of the interface. In agreement with recent theoretical investigations (Eur. Phys. J. E 18 407 
(2005)) we find a strongly anistropic capillary wave spectrum. For the wave-numbers accessed in 
our simulations, the spectrum is quadratic, i. e. elasticity does not play a role. We conjecture that 
this effect is due to the strong bending rigidity of the director field in suspensions of spherocylinders. 

PACS numbers: 61.20. Ja,64.70.Md,64.70.Ja 



I. INTRODUCTION 

Anisotropic particles form liquid crystals at low enough 
temperatures or high enough pressures. Elongated par- 
ticles, in particular, undergo a phase transition from a 
phase in which particle orientations and positions are dis- 
ordered (isotropic phase) to a phase in which orientations 
are aligned and positions are disordered (nematic phase) . 
In many materials the direction of preferred alignment of 
the particles (the director n) can be easily manipulated 
with electric or magnetic fields, which means that the 
optical properties of these materials can be tuned. This 
has made liquid crystals the basis for a large range of 
technological devicesi. 

Liquid crystals also pose fundamental questions about 
the role which particle anisotropy plays in phase tran- 
sitions and interfacial properties. In the 1940s Onsager 
showed that the transition between the isotropic phase (I) 
and the nematic phase (N) is of entropic nature and that 
it can be explained by a simple geometrical argument^. 
However, the calculation of interfacial properties in the 
Onsager model proved to be far more difficult than the 
prediction of the bulk phase transition. 

A. Recent work on mean field properties of the IN 
interface 

Recently there have been several theoretical investiga- 
tions on the IN-interfaco^i'^i^i'^i'^i^i^d°i". In the frame- 
work of density functional theory (DFT) all these stud- 
ies showed that the interfacial tension of hard rods is 
minimal in the case of parallel anchoring, i. e. when the 
nematic director is parallel to the plane of the inter- 
face. Also, alignment is induced in the isotropic phase 
close to the interface. McMuUen^ showed by means of 
a variational method based on a square gradient expan- 



sion for the excluded volume that the density profile is 
monotonic in the case of parallel alignment. This result 
was confirmed by the full numerical studies of Chen and 
Noolandi^'^. They also found small biaxial effects at the 
interface, and a value of ji^LD/ ksT — 0.187 ± 0.001 
for the interfacial tension 71N, where L is the length of 
the particle and D its diameter. This value was 50% 
lower than the result of the variational study by Mc- 
MuUen. Shundyak and van Roij^fi developed a numeri- 
cal method for a full description of the excluded volume 
based on the assumption that the biaxiality is small, as 
Chen and Noolandi previously suggested. This study 
showed indeed a small biaxiality, being thus consistent 
with its hypothesis. The effect of using a finer spa- 
tial grid was that an even smaller interfacial tension 
value j^LD/ksT = 0.156 ± 0.001 was found (which 
agrees with recent simulational data obtained by Vink 
and coworkersi^ii^) . The density and nematic order pa- 
rameter profiles were monotonic as Chen and Noolandi 
had also observed. Only a small non-monotonicity ap- 
peared in the density at the isotropic side in the case 
of normal anchoring, i.e. when the nematic director is 
normal to the plane of the interfac o^i^d° . 

Allen and coworkers have performed several simula- 
tion studies of closely related systems: they compared 
Onsager DFT with simulations of soft rodsM and with 
simulations of hard ellipsoids confined between walls—; 
further they computed the interfacial tension of the IN- 
interface in hard and soft ellipsoids from the pressure 
tensor profile^®, and studied capillary waves on the IN- 
interface of repulsive ellipsoidsii. In these studies they 
found that the interfacial tension agrees well with On- 
sager DFT, while the absolute coexistence densities are 
overestimated by theory. 

In another simulation study nematic wetting films be- 
tween a hard wall and an isotropic fluid of hard rods 
have been observedi^. So far, however, the recent mean 
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field predictions on the structure and tension of the IN- 
interface in long cylindrical rods have not been tested 
against computer simulations. 

If one sets out to compare mean-field results to com- 
puter simulations, the role of fluctuations needs to be 
considered. In the following section we briefly summa- 
rize the concept of capillary waves, and in section Fl CI we 
discuss how we compare mean-field "intrinsic" profiles to 
our simulation results. 



B. Capillary waves 

At finite temperature, an interface is never planar be- 
cause of thermal fluctuations ("capillary waves"). As 
each deviation from the planar shape costs energy, the 
capillary wave spectrum is determined by the balance be- 
tween thermal energy and interfacial tension (and other 
contributions to the free energy cost of the interface, such 
as bending rigidity.) In the simplest case, the capillary 
wave spectrum is goverend by the free energy cost of the 
interfacial area, which is added due to the undulations. 
Assuming that the interface position can be parametrized 
by a function h{x,y), and that the local distortions are 
small, the wave spectrum is predicted to be 

(IMq)l^> = ^ , (1) 

where ft.(q) is the Fourier transform of the local position 
h(x, y) of the fluctuating surface, q is the wavevector and 
7 is the interfacial tension. Note that this expression can- 
not hold for (? — > 0, because the amplitude would diverge. 
In real systems, gravity or finite system size, for example, 
introduce a cutoff at small q. For an introduction to the 
subject of capillary waves see for example the book by 
Safrani^. 

However, the capillary wave spectrum of the IN- 
interface cannot be described in this simple manner. 
Liquid crystals exhibit long-range elastic interactions, 
which keep the IN interface anisotropic up to very large 
length scales. (The interface becomes isotropic only on 
scales larger than the length over which deformation of 
the director field costs less energy than k-QT.) Also, the 
interfacial tension 71N is not a constant, but a function 
of the angle (j) between the director and the interface. 
Therefore, fluctuations parallel to the director will pro- 
duce different free energy costs than fluctuations perpen- 
dicular to the director. This separates the IN-interface 
from most other interfaces for which capillary waves have 
been studied. 

Recently there have been several theoretical studies 
on waves on the IN-interface22iSiiSS. In particular, the 
IN-capillary wave spectrum has been analyzed by Elgeti 
and Schmid within Frank elastic theory and Landau-de 
Gennes theorjiSS. For the Landau-de Gennes approach 
their study shows that 



where 7]^ is the interfacial tension for in-plane alignment 
(i. e. (/i = 0), qy is the component of the unit vector q/q 
in the direction of the bulk director and 6 is a constant, 
which contains information on the elasticity of the direc- 
tor field and the anisotropy of the interfacial tension. (A 
more detailed expression can be found in the article by 
Elgeti and Schmidt.) Thus, the spectrum is isotropic in 
the leading g^-term, but shows a strong anisotropy in the 
next order g^-term. An anisotropic contribution to the 
q^-term, which one would expect to find because of the 
anisotropic interfacial tension, is removed by the elastic 
interactions. Both the theoretical work and the computer 
simulation on repulsive ellipsoids by Akino, Schmid and 
Allenii show that the amplitudes of the waves are largest 
in the direction perpendicular to the director. Here, we 
compare these results to the hard-spherocylinder case. 



C. Comparing DFT and simulation 

The IN-transition in rods is of first order and it in- 
creases in strength with increasing aspect ratio. For large 
aspect ratios the only relevant fluctuations in a finite 
simulation-box are therefore capillary waves. For the 
further discussion, we assume that there is an "intrin- 
sic" mean-field interface, which is broadened by capillary 
waves in a manner that is statistically uncorrelated to 
the properties of the "intrinsic" interface. This concept 
has been under debate, because it clearly fails for short 
wavelengths i. e. wavelengths of the order of the size 
of the particlesSiSiiSSiSSiSl. However, it has proven to 
be successful in studies for which the short wavelength 
limit was not relevant (see for example22i22i2£) . Here, we 
are interested in a comparison of simulations to Onsager 
theory and Landau-de Gennes theory, i. e. we are not 
studying molecular details and therefore we assume the 
intrinsic profile concept to be valid. 

We test which effect fluctuations have on our DFT re- 
sults, by convoluting the "intrinsic" DFT profile fp^"^ [z) 
with a simple Ansatz for the capillary wave spectrum 

2s2 
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where h is the local height of the interface, is the local 
mean-square displacement of the interfacial height 

s — - — m 

27r7 V 9min 

Qmax is the length of the maximum wavevector (given 
by the particle diamter) and gmin is the length of the 
minimum wavevector (given by the size of the box) . The 
estimated apparent profile then becomes 



dh p^^'^ {z - h)P{h) 



IN 



(IMq)l= 



[bq^ + Oiq" 



This Ansatz for the capillary wave spectrum certainly 
does not capture the specific properties of the IN- 
interface described in the previous section. But we will 
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show in the results section, that it suffices to test whether 
a quantity can be compared between DFT and simulation 
or not. 

The article is structured as follows: First we introduce 
the model and the simulation techniques. Then we dis- 
cuss the interfacial profiles and compare them to DFT 
results. In the third part we present the analysis of the 
capillary wave spectrum and in the last section we sum- 
marize. 



II. MODEL AND ORDER PARAMETERS 
A. DFT 

The basic ingredient of density functional theory is the 
grand potential i^[p\ given as a functional of the one- 
particle distribution p{r, u), with r the center-of-mass of 
a rod and u the unit vector of the long axis of a rod: 

m = J drdup(r, u) (log[p(r, u)L^D] - 1 - /3m) 

~ \j drdudr'du7(r,u;r',u')/3(r,u)p(r',u'), 

where f3 = l/fc^T, and / is the Mayer function which 
equals — 1 if the rods overlap and otherwise^. Here we 
focus on long hard spherocylinders with a length L much 
larger than the diameter D. In the limit D/L — > 0, in 
which we perform the DFT calculations, Onsager's sec- 
ond virial functional is expected to be accurateSiS*i2*ii. 
For a given chemical potential /x the equilibrium den- 
sity minimizes the functional, and the minimum condi- 
tion 6il[p]/Sp{r,u) = gives rise to a nonlinear integral 
equation for the equilibrium profile p(r, u) : 

log[p{r,u)L^D]~ J dr'du7(r,u;r',u')p(r',u') 

The minimum value of the functional, which is obtained 
after insertion of the equilibrium profile into the func- 
tional, gives the equilibrium grand potential —pV for a 
bulk system of volume V at pressure p, and —pV + jA 
for a system with a planar interfacial area A (with 7 the 
interfacial tension). In the latter geometry we assume 
that the distribution is independent of the in-plane coor- 
dinates X and y, such that the Euler-Lagrange equation 
is to be solved for p{z,0,(fi), with z the cartesian coor- 
dinate perpendicular to the interface and 6 and ip the 
polar and azimuthal angle of u with respect to the ne- 
matic director far from the interface, respectively. Here 
we impose boundary conditions such that the coexist- 
ing isotropic and ncmatic bulk phases are obtained at 
\z\ — > 00. From the equilibrium profile one extracts^i^iifi 
the total density and the order parameter profile as 
p{z) = J du p {z, u) and S{z) ~ J du P2{n-u)p{z,u)/ p{z), 
respectively, with P2{x) = (3x^ — l)/2 the second Leg- 
endre polynomial. The biaxiality profile is given by 
a{z) = (3/2) J du p{z , u) siii^ 9 cos{2ip) / p{z) . We present 



the numerical results of the corresponding profiles in 
Figs. niEl andO An accurate numerical methodic gives 
for the surface tension -y^'"'^ LD /ksT = 0.156 ± 0.001 at 

parallel anchoring, and -i^^'^LD/ksT = 0.265 ± 0.001 
at normal anchoring. 

B. Simulations 

As in the Onsager model, we consider spherocylinders, 
each consisting of a cylinder of length L and diameter D 
capped by two hemispheres of diameter D. The inter- 
action between two particles i and j is given by a pair 
potential of the form 

v -ir) ~ 

I otherwise, 

where r is the distance between the cylinder axes. 

To investigate the IN transition, the density and the 
average rod alignment are used as order parameters. 
We denote the density as the dimensionless quantity 
p* = [1: /A) L'^D N /V where N is the number of parti- 
cles, and V is the volume of the simulation box. The 
average alignment is defined in terms of the orientational 
traceless tensor Q with the elements 

1 ^ 

Qafi = ^ ^ {^UiaUip - Safl) , 
i=l 

where is the a component (a = x, y, z) of the unit 
vector along the axis of particle i and Sap the Kronecker 
delta. Diagonalization yields three eigenvalues, 

S, -{S-a)/2, -(5 + a)/2, 

where S is the nematic order parameter and a is the 
biaxiality, which is equivalent to the definition used in 
the DFT calculations. 



III. SIMULATION METHOD 

Simulations were performed in the canonical ensemble 
using an elongated box with periodic boundary condi- 
tions and dimensions x Ly x = 10 x 10 x 20 L^. 
By rotational and translational Monte-Carlo moves phase 
space was accessed. The aspect ratio of the rods was 

= 15. The choice of L/D=15 is motivated by the fol- 
lowing requirements: on the one hand, the IN transition 
becomes weaker with decreasing L/D. At short L/D fluc- 
tuations in the bulk phases come into play, making the 
comparison with DFT very difficult. On the other hand 
the CPU time required to detect an overlap between two 
rods increases with L/D. Thus a compromise has to be 
made. In a previous study on ellipsoids, L/D=15 had 
been chosenii. We used the same to be able to compare 
our results to this study. 
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The system was prepared in isotropic-nematic coex- 
istence and the interfacial planes were located in the 
(x, ?/)-plane. Due to the periodic boundary conditions 
at least two interfaces were required. The overall density 
was chosen between the coexistence values, which we had 
previously determined by grand canonical simulationai^. 
With this geometry about 10^ particles were required. 
Parallel interfacial anchoring as well as normal anchor- 
ing was investigated. 

Normal anchoring is metastable, therefore the direc- 
tor would rotate towards a parallel orientation during a 
sufficiently long simulation. The time-scale necessary for 
such a rotation is given by the time it takes to cross the 
free-energy barrier for the creation of a defect contain- 
ing a reoriented region. As the system is relatively stiff 
with respect to director field bending, we expect this bar- 
rier to be rather high. During our simulation the mean 
z-component of the director decreased by less than half 
a percent, while the fluctuations around this value were 
much larger than that. We have produced some data sets 
using a constraint on the orientation of the director, as it 
was introduced by Allen and coworkers"^^. However, the 
results did not differ from the ones obtained in uncon- 
straint simulations. 

After the equilibration run of about 10'' MC sweeps 
(i. e. moves per particle) the values of the density and the 
orientational order parameter far away from the interface 
did not change anymore. The bulk values were reached 
in clear plateaus and they were in agreement with the 
results from the grand canonical simulations. These con- 
ditions were used as evidence for an equilibrated state. 
For data production about 1.5 x 10^ MC sweeps were 
performed. 
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FIG. 1: Density and order parameter profiles for parallel an- 
choring. 

Upper graph: Simulation results of the reduced density p* 
(triangles) and the nematic order parameter S (squares). 
Dashed and solid lines are tanh-fits, which serve to guide the 
eye. The bulk values of the density according to simulations 
in the grand-canonical ensemble— are indicated by a dotted 
line. The shift S between the inflection points of the profiles 
comes out clearly. Inset: Decay of the density and the order 
parameter to their isotropic bulk values on a log scale. Both 
quantities have the same decay length. 

Lower graph: Results from DFT calculations for density 
(dashed line) and nematic order prameter (solid line) pro- 
files. 

Note that Onsager DFT overestimates the coexistence den- 
sities, while the shift S between the interfaces agrees well. 
Differences between DFT and simulations due to fluctuations 
are discussed in the text. 



IV. RESULTS 

A. Interfacial profiles 

To investigate the interfacial profiles, the box was de- 
vided into 60 slices in z— direction, and the density, the 
nematic order parameter, and the biaxiality order pa- 
rameter were sampled in each slice. By averaging over 
3 X 10^ profiles the statistical error on the density, which 
is much larger than the error on the order-parameter, be- 
came reasonably small and a clear biaxiality signal could 
be obtained. 

Note, that the total average volume of nematic phase 
is conserved, but that the center of the nematic layer can 
shift along the z— direction without creating any ener- 
getic or entropic costs. To account for the freedom of 
the system to shift in the z— direction, the profiles had 
to be centered before sampling. We shifted them such 
that their centers of mass lay in the middle of the simu- 
lation box. (Fluctuations of the total amount of nematic 
phase lead to broadening of the averaged profiles^Si^. 
This effect will be discussed in section llVBI 'l The result- 
ing profiles for parallel anchoring from theory as well as 
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FIG. 2: Biaxiality profile from simulations (solid line) and 
DFT calculations (dashed line) for parallel anchoring. 



simulation are plotted in Figs, ^a-ndl^l 

Figs. Hand 121 show that the interface induces both ne- 
matic order and biaxiality on the isotropic side near the 
interface. This is in agreement with DFT calculations. 
The absolute values of the densities, however, differ con- 
siderably. This effect has been observed in similar sys- 
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FIG. 3: Density and order parameter profiles for normal an- 
clioring. 

Upper grapti: Reduced density p* (triangles) and nematic 
order parameter S (squares). Dashed and solid lines are 
tanh-fits, which serve to guide the eye. The bulk values of 
the density according to simulations in the grand-canonical 
ensembleW^ are indicated by a dotted line. 
Lower graph: DFT calculations for density (dashed line) and 
order parameter (solid line) profiles. Both in simulation and 
theory, the density profile shows a small non-monoticity on 
the isotropic side. 



terns by Allen and coworkersiai^. The induced nematic 
order on the isotropic side is reflected in the shift be- 
tween the inflection points of the density profile and the 
order parameter profile, (5 = 0.38(0) L. DFT gives a sim- 
ilar result, 6 — 0.451"^"'. Also in computer simulations 
of soft spherocylinders 8 = 0.37(4) L was found^'*. For 
macroscopic system sizes, this quantity is not a constant 
because it is affected by capillary waves. However, the 
profiles are almost symmetric with respect to their in- 
flection points. Therefore convolution with a symmetric 
capillary wave spectrum Eq.(|2I) will hardly shift the in- 
flection points and thus will not lead to changes in b 
unless very large systems are considered. We tested this 
both, by changing the size of the simulation box and by 
computing the apparent profiles p'^^^ij.) and S'°^*(z) from 
the DFT profiles. The changes in b were of the order of 
10%. 

DFT also predicts that the density and order param- 
eter profiles should decay on the same length-scale, be- 
cause the correlation length ^ is the only characteristic 
length scale for these decays. The inset of Fig. shows 
simulation results for both profiles on a logarithmic scale 
revealing agreement with the theoretical prediction. (For 
the system size used in the simulation, this property is 
not affected by capillary waves, because it is measured 
far away from the interface.) 

Fig. 13 shows the profiles for normal anchoring. Again 
there is good agreement between simulation and theory. 



In particular, both show a weakly non-monotonic be- 
haviour of the density profile on the isotropic side. For 
increasing system-sizes this dip is smoothened out by in- 
terfacial broadening. 

The present simulations confirm the theoretical 
claims^iiSiii that, for the homeotropic (biaxial) IN- 
interface, the density and order parameter profiles are 
monotonic and that the interfacial biaxiality is small. 
For the uniaxially symmetric IN-interface the simulations 
confirm the theoretically predicted small non-monotonic 
feature in the density profile. 



B. Blocking analysis 

In order to analyze capillary waves the interfacial po- 
sition has to be extracted from the simulation data. This 
requires the choice of a definition of what constitutes the 
interface. Different choices of definition can lead to differ- 
ences in the capillary wave spectra, in particular for short 
wavelengths. The simplest approach is to cut the system 
into blocks of length and cross section area B x B 
and to compute the Gibbs dividing surface for each block. 
However, the smaller the block, the larger are the fluctua- 
tions of the densities of the coexisting phases around their 
bulk values. This problem is circumvented by the method 
we describe below and similar methodsiiffi which are 
based on sampling the differences between local profiles 
for each block in a configuration and the overall profile 
of the same configuration. These methods do not take 
into account the molecular details of the interface and are 
therefore not useful for studies of the capillary wave spec- 
trum at short length-scales. If one is interested in these, 
the pivot particle method by Tarazona and Chacon^ can 
be used, which proceeds via a different route: First par- 
ticles, which form the outmost layer of one phase, are 
identified. Then the minimal surface through the posi- 
tions of these particles is constructed. Iteratively, more 
particles from the same phase which lie close to this sur- 
face are added and a new minimal surface is constructed 
until all surface particles are incorporated. This method 
has shown to produce considerable differences in capil- 
lary wave spectra at short wavelengths, and it would be 
worthwhile applying it to the IN-interface. However, as 
we are interested in a comparison with Landau-de Gennes 
theory, a continuum theory, we use a method that is less 
accurate on the molecular level. 

In order to analyze the capillary waves, we used the fol- 
lowing procedure: We first determined the average plane 
of the interface (averaged over the whole box cross sec- 
tion) and then measured undulations around this refer- 
ence plane. It is, however, not sensible to determine a 
reference plane for each individual system snapshot, as 
the position of the averaged profile is itself a fluctuating 
quantity. Instead we used averages over 10^ MC sweeps 
to determine stable reference planes. For each snapshot, 
we shifted the system such that its center of mass was 
in the center of the box. Then the box was divided into 
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FIG. 4: Sketch of the blocking scheme. Details can be found 
in the text. 
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FIG. 5: Interfacial profiles for various block sizes. The inset 
shows a close-up view with a logarithmic scale of the profile 
near the center. 



two parts, each containing one of the two interfaces. The 
interfaces were then investigated separately. In order to 
determine the reference plane of each interface we fitted 
the data with 



S(z) 



Si 



Si) tanh 



zq 



ws 



(3) 



where Sn and Si are the bulk nematic and isotropic order 
parameters respectively. The width ws and the position 
of the point of inflection zq are fit parameters. Obvi- 
ously, the profiles will not be perfect hyperbolic tangent 
functions, because of the asymmetry in orientational or- 
der between the two phases. However, the fit only serves 
to determine an average interfacial position as a basis for 
the following analysis: As sketched in Fig. 01 the box was 
then devided into square blocks in the (x, ?/)-plane with 
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FIG. 6: Squared interfacial width versus ln(i3) for parallel 
anchoring (solid line) and normal anchoring (dotted line). 



side lengths B varing from L.j;/7 to L^- (Note that even 
the smallest block has a side length, which is larger than a 
particle length. Thus we are far away from the molecular 
diameter D.) Each block was then divided into 25 bins in 
the interval [zo — 4L, zq + 4L] in z-direction. The pro- 
files which we measured in these subsystems were again 
fitted by the tanh- function Eq.Q in order to find their 
points of inflection. All subsystem profiles for a given 
block size B were centerend at their points of inflection. 
Then Q was averaged over all profiles and the average 
profiles were fitted by Eq.Q again in order to sample 
the squared interfacial width w| and the interfacial posi- 
tions z(xo,?/o) (where XQ,yo denotes the position of the 
block center). 255 system snapshots, i. e. 510 snapshots 
of interfaces, have been used for the blocking analysis. 
The smallest subsystem contained about 7500 particles, 
which was sufficient to find a clear order parameter sig- 
nal. 

In order to make sure that the snapshots were decor- 
related, we checked the autocorrelation functions of Wg 
and z for various block sizes. We also ruled out possible 
effects due to the finite bin size by comparing results for 
different discretisations. The size of a bin is determined 
by the following compromise: On the one hand a small 
bin size provides a good resolution of the profile, but on 
the other hand the computation of order parameter aver- 
age requires at least 100 particles. For bin sizes between 
10 D and 20 D the interfacial width as well as the com- 
puted interfacial tension remained constant within the 
numerical accuracy. 

Fig. shows the averaged order parameter profiles 
for various block sizes. If a capillary wave spectrum is 
of the form of Eq.(^) the interfacial width behaves as 
Ws oc y/\n{B)^^. Fig. shows Wg vs. \n{B) for paral- 
lel and normal anchoring. For parallel anchoring we find 
the predicted linear dependence. For normal anchoring 
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FIG. 7: Inverse mean-squared fourier components of the in- 
terface position versus wave number squared (m = q{L^/2n)) 
for parallel anchoring. 
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FIG. 9: L^/ ((|/i(q)|^)m^). For m > 2 amplitudes are pro- 
portional to The anisotropy of the interfacial tension 
shows already on the q^-level. 
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FIG. 8: Inverse mean-squared fourier components of the in- 
terface position versus wave number squared for normal an- 
choring. 



there might be deviations. However, the accuracy of the 
simulations does not allow for a more detailed statement. 
The error-bars have been determined by the "jackknife- 
method"^. 



C. Analysis of the capillary wave spectrum 

A similar analysis as described in the previous section 
was used in order to measure the capillary wave spec- 
trum. For blocks of size x Ly, we sampled the 
distances between the points of inflection of the local pro- 
files and the point of inflection of the reference profile. In 
order to test the stability of our results, we also applied 
the method described by Akino and coworkers^^. The ac- 



curacy was less good, but the results remained the same. 
Further we checked that the distribution of interfacial 
heights P{h) was Gaussian, indicating that the system 
was properly thermalized. Correlations between the un- 
dulations decayed within 10* MC sweeps for the larger 
values of q and 5-10'' for the smallest value. Hence we 
are confident that the snapshots over which we averaged 
were sufficiently decorrelated. 

Figs. and |S1 show L^/(|/i(q)p) vs. wave num- 
ber squared for parallel and normal anchoring. (We 
have plotted m = q(Lx/27r) instead of q to facili- 
tate direct comparison to the results of Akino and co- 
authors^''.) Again the error-bars have been determined 
by the "jackknife-method" . In the case of parallel an- 
choring the director was oriented along the y-axis. As ex- 
pected, the spectrum is anisotropic for parallel anchoring 
and isotropic for normal anchoring. For parallel anchor- 
ing, we find in agreement with theorjiSS that fiuctuations 
parallel to the director (dotted line) are smaller than fiuc- 
tuations perpendicular to the director (solid line). 

The magnitude of the fluctuations compares well with 
the results Akino and co-authors have obtained for 
ellipsoidsii. However, their results deviate from the q^ 
dependence already for ^ 10. To check whether there 
are contributions of higher order than g^, we have plotted 
L^/ ((|ft,(q)p)TO^) versus wave number for parallel anchor- 
ing in Fig. El For m > 2, both components are almost 
constant. (The large values for very small wave numbers 
are due to the finite system size.) Parallel to the direc- 
tor, there might be small deviations from the q^ depen- 
dence, but the results are not accurate enough to allow 
for detailed conclusions. Additionally we find that there 
is a difference in the absolute value L^/ ((|/i(q)p)TO^) 
between parallel and perpendicular contributions. This 
means that for the regime of wave-numbers studied in 
this work, the anisotropy of the interfacial tension does 
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FIG. 10; Sketch of capillary waves: 

Left: Strong interfacial anchoring — > anisotropy of the spec- 
trum is dominated by elasticity 

Right: Strong bending rigidity — > anisotropy of the spectrum 
is dominated by the anisotropy of the interfacial tension 



enter already on the g^-level. 

Fig-Elshows a sketch of capillary waves in two extreme 
cases: in the case on the left, parallel anchoring is so 
strong that the particles' orientations have to follow the 
shape of the interface. The anisotropy of the spectrum is 
then dominated by contributions from bending elasticity. 
In the case on the right, the bending rigidity is so strong, 
that the anisotropy of the spectrum is dominated by con- 
tributions from the anisotropy of the interfacial tension. 
From the fact that we do not observe contributions due 
to bending, but strong contributions due to interfacial 
tension anisotropy, we conclude that spherocylinders fall 
into the category sketched on the right. 

V. DISCUSSION AND SUMMARY 

We have presented a comparative study of computer 
simulations and Onsager DFT for the isotropic-nematic 



interface in hard spherocylinders. Although Onsager the- 
ory overestimates the isotropic-nematic coexistence den- 
sities, interfacial properties such as the the shift between 
the density profile and the nematic order parameter pro- 
file are in agreement. This observation accords with the 
studies by Allen and coworkers on ellipsoidai^ and soft 
spherocylindersi^. 

Our simulations confirm that, for the case of in-plane 
orientation of the director at the interface, the density 
and order parameter profiles are monotonic and that the 
interfacial biaxiality is small. For the case of normal ori- 
entation the simulations agree with the theoretically pre- 
dicted small non-monotonic feature in the density profile. 

The capillary wave spectrum of the interface (in 
the case of parallel anchoring) is anisotropic, as ex- 
pected. In contrast to the case of ellipsoids, we find that 
L^/(|ft,(q)P) (X for the entire range of wave vectors 
which could be accessed in this study. Also we observe 
anisotropic contributions already on the g^-level. We in- 
terpet both facts as signals of the a very stiff director 
field, which prohibits bending along the fluctuating in- 
terface. 



Acknowledgments 

We thank the Deutsche Forschungsgemeinschaft 
(DFG) for support (TR6/D5) and K. Binder, R. Vink, 
P. van der School and J. Horbach for helpful suggestions. 
TS is supported by the Emmy Noether program of the 
DFG and by the MWFZ. Allocation of computer time 
on the JUMP at the Forschungszentrum Jiilich is grate- 
fully acknowledged. This work is part of the research 
program of the "Stichting voor Fundamcntccl Onderzoek 
der Materie (FOM)" , which is financially supported by 
the "Nederlandse organisatie voor Wetenschappelijk On- 
derzoek (NWO)". 



^ P. G. de Gennes and J. Prost, The Physics of Liquid Crys- 
tals (Clarendon Press, Oxford, 1993). 

2 L. Onsager, Ann. N. Y. Acad. Sci. 51, 627 (1949). 

^ W. E. McMuUen, Phys. Rev. A 38, 6384 (1988). 

" Z. Y. Chen and J. Noolandi, Phys. Rev. A 45, 2389 (1992). 

^ Z. Y. Chen, Phys. Rev. E 47, 3765 (1993). 

" D. L. Koch and O. G. Harlen, Macromolecules 32, 219 
(1999). 

'' P. van der School, J. Phys. Chem. B 103, 8804 (1999). 
* E. Velasco, L. Mederos, and D. E. Sullivan, Phys. Rev. E 

66, 021708 (2002). 
^ M. P. Allen, Chem. Phys. Lett. 331, 513 (2000). 

K. Shundyak and R. van Roij, J. Phys.: Condens. Matter 

13, 4789 (2001). 
" K. Shundyak and R. van Roij, Phys. Rev. E 68, 061703 



(2003). 

R. L. C. Vink and T. Schilling, Phys. Rev. E 71, 051716 
(2005). 

R. L. C. Vink, S. Wolfsheimer, and T. SchiUing (2005). 
M. S. Al-Barwani and M. P. Allen, Phys. Rev. E 62, 6706 
(2000). 

M. P. Allen, J. Chem. Phys. 112, 5447 (2000). 

A. J. McDonald, M. P. Allen, and F. Schmid, Phys. Rev. 

E 63, 010701(R) (2000). 

N. Akino, F. Schmid, and M. P. Allen, Phys. Rev. E 63, 
041706 (2001). 

M. Dijkstra, R. van Roij, and R. Evans, Phys. Rev. E. 63, 
051703 (2001); R. van Roij, M. Dijkstra, and R. Evans, 
Europhys. Lett. 49, 350 (2000); M. Dijkstra, and R. van 
Roij, J. Phys. Cond. Matt. 17, S3507 (2005). 



9 



S. A. Safran, Statistical Thermodynamics of Surfaces, 
Interfaces, and Membranes (Addison- Wesley Publishing 
Company, New York, 1994). 

V. Popa-Nita and T. J. Sluckin, Phys. Rev. E 66, 041703 

(2002) . 

V. Popa-Nita and P. Oswald, Phys. Rev. E 68, 061707 

(2003) . 

J. Elgeti and F. Schmid, Eur. Phys. J. E 18, 407 (2005). 
P. Tarazona and E. Chacon, Phys. Rev. B 70, 235407 

(2004) . 

J. Stecki, J. Chem. Phys. 109, 5002 (1998). 

M. Napiorkoswski and S. Dietrich, Phys. Rev. E 47, 1836 

(1993). 

K. Mecke and S. Dietrich, Phys. Rev. E 59, 6766 (1999). 
A. Robledo, C. Vaxea, and V. Romero-Ronchin, Physica A 
177, 474 (1991). 



® A. Werner, F. Schmid, M. Miiller, and K. .Binder, J. Chem. 

Phys. 107, 8175 (1997). 
'^^ A. Werner, F. Schmid, M. Miiller, and K. .Binder, Phys. 

Rev. E 59, 728 (1999). 
^° R. Vink, J. Horbach, and K. Binder, J. Chem. Phys. 122, 

134905 (2005). 

M. P. Allen, M. A. Warren, M. R. Wilson, A. Sauron, and 

W. Smith, J. Chem. Phys. 105, 2850 (1996). 

H. L. Topper and W. J. Briols, J. Chem. Phys. 116, 5186 

(2002). 

K. Binder, Phys. Rev. A 25, 1699 (1982). 

K. Shundyak, Dissertation, University of Utrecht, The 

Netherlands (2004). 
'^'^ B. Efron, Annals of Statistics 7, 1 (1979). 



